Haciendo un análisis de la serie apartir de 1948 hasta 2020 Marzo. A primera vista se observa un crecimeinto en el año 2008 causado por la crisis del 2008, por lo que se observa un crecimiento en la tasa de desempleo. Crisis que se apartir del 2010 comienza su descenso. De igual manera se observa un
## Rows: 874
## Columns: 2
## $ time <dttm> 1948-01-01, 1948-02-01, 1948-03-01, 1948-04-01, 1948-05-0...
## $ UnemRate <dbl> 4.0, 4.7, 4.5, 4.0, 3.4, 3.9, 3.9, 3.6, 3.4, 2.9, 3.3, 3.6...
Descomponemos la serie para observar su componentes: tendencia ,estacionalidad y los residuales
Se observa que la varianza es homogendea por l oque le aplicamos una trnadofrmación. Ocupando el metodo de guerrero obtenemos una Lambda de .99.
## Registered S3 method overwritten by 'quantmod':
## method from
## as.zoo.data.frame zoo
## -- Attaching packages ---------------------------------------------------------------------------------- fpp2 2.4 --
## v forecast 8.13 v expsmooth 2.3
## v fma 2.4
##
Gráficando las ACF Y PACF se observa que necesitamos una diferencia, ya que se presenta una persistencia muy alta.
APlicando uan diferencia
Se observa un componente estacional a cada doce periodos lo que nos indica que podria ser un Sarima(1,1,1)(1,0,0)12. Ahora, cada tercer rezago se observa un incremento, continuado de un descenso, por lo que en la parte de estacionariedad propondriamos un modelo ARIMA(3,1,2) o (3,1,0) De esta manera proponemos dos modelos: SARIMA(3,1,2)X(1,0,0)12 Y SARIMA(3,1,0)X(1,0,0)12
Estimamos dos modelos SARIMA:
El primero SARIMA(3,1,2)X(1,0,0)12
##
## Call:
## arima(x = ts_des, order = c(3, 1, 2), seasonal = list(order = c(1, 0, 0), period = 12))
##
## Coefficients:
## ar1 ar2 ar3 ma1 ma2 sar1
## -0.8126 0.3399 0.2669 -0.0487 -0.9513 0.8386
## s.e. 0.0352 0.0422 0.0340 0.0145 0.0145 0.0200
##
## sigma^2 estimated as 0.07592: log likelihood = -120.49, aic = 254.98
##
## Training set error measures:
## ME RMSE MAE MPE MAPE MASE ACF1
## Training set -0.001657906 0.275373 0.2084416 NaN Inf 0.4437551 -0.008711197
El segundo modelo observado: SARIMA(3,1,0)X(1,0,0)12
mod2 <- arima(x = ts_des, order = c(3, 1, 0), seasonal = list(order = c(1, 0, 0), period = 12))
summary(mod2)
##
## Call:
## arima(x = ts_des, order = c(3, 1, 0), seasonal = list(order = c(1, 0, 0), period = 12))
##
## Coefficients:
## ar1 ar2 ar3 sar1
## -0.7514 -0.3999 -0.1795 0.8445
## s.e. 0.0340 0.0402 0.0338 0.0182
##
## sigma^2 estimated as 0.08786: log likelihood = -183.32, aic = 376.63
##
## Training set error measures:
## ME RMSE MAE MPE MAPE MASE ACF1
## Training set 0.0004268723 0.2962325 0.2240249 NaN Inf 0.4769307 -0.03826402
estos dos modelos el que tiene el menor grado de AIC es el primero.
Observando la prueba Diebold-Mariano se observa que el p-value de de la última prueba es mayor a cero por lo qe se rechaza la hipotesis nula, donde el estadistico de prueba es diferente de cero. Por lo que se acepta la hipotesis alternativa de que el mod1 es
library(forecast)
dm.test(residuals(mod1),residuals(mod2),alternative = c("two.sided"),h=1)
##
## Diebold-Mariano Test
##
## data: residuals(mod1)residuals(mod2)
## DM = -3.9309, Forecast horizon = 1, Loss function power = 2, p-value =
## 9.137e-05
## alternative hypothesis: two.sided
dm.test(residuals(mod1),residuals(mod2),alternative = c("less"),h=1)
##
## Diebold-Mariano Test
##
## data: residuals(mod1)residuals(mod2)
## DM = -3.9309, Forecast horizon = 1, Loss function power = 2, p-value =
## 4.569e-05
## alternative hypothesis: less
dm.test(residuals(mod1),residuals(mod2),alternative = c("greater"),h=1)
##
## Diebold-Mariano Test
##
## data: residuals(mod1)residuals(mod2)
## DM = -3.9309, Forecast horizon = 1, Loss function power = 2, p-value =
## 1
## alternative hypothesis: greater
Graficaremos los Pronosticos para observar de los deos modelos
mod1_res<- sarima.for(ts_des, 7,3,1,2, 1,0,0,12)
mod1_res<- sarima.for(ts_des, 7,3,1,0, 1,0,0,12)
Se observa que -22.59 esta en el area de Rechazo por lo que suponemos que no existe Raíz unitaria. Observando el estadistico tt se observa que no es significativo de cero. por lo que hacemos lel modelo sin tendecia lineal.
library("urca")
ad.mod1 <- ur.df(ts_des, type="trend", selectlags="AIC")
summary(ad.mod1)
##
## ###############################################
## # Augmented Dickey-Fuller Test Unit Root Test #
## ###############################################
##
## Test regression trend
##
##
## Call:
## lm(formula = z.diff ~ z.lag.1 + 1 + tt + z.diff.lag)
##
## Residuals:
## Min 1Q Median 3Q Max
## -1.1606 -0.2690 -0.1023 0.1814 2.3118
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 9.178e-03 3.096e-02 0.296 0.766971
## z.lag.1 -1.035e+00 4.583e-02 -22.593 < 2e-16 ***
## tt -2.095e-05 6.191e-05 -0.338 0.735163
## z.diff.lag 1.290e-01 3.382e-02 3.816 0.000145 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.4538 on 860 degrees of freedom
## Multiple R-squared: 0.4669, Adjusted R-squared: 0.4651
## F-statistic: 251.1 on 3 and 860 DF, p-value: < 2.2e-16
##
##
## Value of test-statistic is: -22.5931 170.1521 255.2273
##
## Critical values for test statistics:
## 1pct 5pct 10pct
## tau3 -3.96 -3.41 -3.12
## phi2 6.09 4.68 4.03
## phi3 8.27 6.25 5.34
#dm.test(residuals(mod1),residuals(mod2),alternative = c("two.sided"),h=7)
Se observa el mismo comportamiento, por lo que volvemos a rechazar que exista la raiz unitaria, y de igual manera 255.43 esta del lado derecho por lo que no es un modelo con deriva.
ad.mod1 <- ur.df(ts_des, type="drift", selectlags="AIC")
summary(ad.mod1)
##
## ###############################################
## # Augmented Dickey-Fuller Test Unit Root Test #
## ###############################################
##
## Test regression drift
##
##
## Call:
## lm(formula = z.diff ~ z.lag.1 + 1 + z.diff.lag)
##
## Residuals:
## Min 1Q Median 3Q Max
## -1.1557 -0.2687 -0.1001 0.1771 2.3140
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 0.0000968 0.0154299 0.006 0.994996
## z.lag.1 -1.0351209 0.0457971 -22.602 < 2e-16 ***
## z.diff.lag 0.1288804 0.0337957 3.814 0.000147 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.4535 on 861 degrees of freedom
## Multiple R-squared: 0.4669, Adjusted R-squared: 0.4656
## F-statistic: 377 on 2 and 861 DF, p-value: < 2.2e-16
##
##
## Value of test-statistic is: -22.6023 255.4336
##
## Critical values for test statistics:
## 1pct 5pct 10pct
## tau2 -3.43 -2.86 -2.57
## phi1 6.43 4.59 3.78
#dm.test(residuals(mod1),residuals(mod2),alternative = c("two.sided"),h=7)
Por último se rechaza con sea un DF sin constante con raíz unitaria, por lo que tenemos un modelo estacionario alrededor de cero
ad.mod1 <- ur.df(ts_des, type="none", selectlags="AIC")
summary(ad.mod1)
##
## ###############################################
## # Augmented Dickey-Fuller Test Unit Root Test #
## ###############################################
##
## Test regression none
##
##
## Call:
## lm(formula = z.diff ~ z.lag.1 - 1 + z.diff.lag)
##
## Residuals:
## Min 1Q Median 3Q Max
## -1.1556 -0.2686 -0.1000 0.1772 2.3140
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## z.lag.1 -1.03512 0.04577 -22.615 < 2e-16 ***
## z.diff.lag 0.12888 0.03378 3.816 0.000145 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.4533 on 862 degrees of freedom
## Multiple R-squared: 0.4669, Adjusted R-squared: 0.4656
## F-statistic: 377.4 on 2 and 862 DF, p-value: < 2.2e-16
##
##
## Value of test-statistic is: -22.6155
##
## Critical values for test statistics:
## 1pct 5pct 10pct
## tau1 -2.58 -1.95 -1.62
#dm.test(residuals(mod1),residuals(mod2),alternative = c("two.sided"),h=7)
Mi modelo se veria muy afectado ante los cambios por el Covid, dado que este comportamiento es nuevo y es algo no predecible. Por lo que mis pronosticos serian incorrectos. Una opción seria tomar pronosticos por Rolling Windows para de esta manera incluir la información de estos últimos meses.